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ABSTRACT 

Synchrotron radiation is commonly observed in connection with shocks of different velocities, rang¬ 
ing from relativistic shocks associated with active galactic nuclei, gamma-ray bursts or microquasars 
to weakly- or non-relativistic flows as those observed in supernova remnants. Recent observations of 
synchrotron emission in protostellar jets are important not only because they extend the range over 
which the acceleration process works, but also because they allow to determine the jet and/or inter¬ 
stellar magnetic field structure, thus giving insights on the jet ejection and collimation mechanisms. 
In this paper, we compute for the first time polarized (synchrotron) and non polarized (thermal X-ray) 
synthetic emission maps from axisymmetrical simulations of magnetized protostellar jets. We consider 
models with different jet velocities and variability, as well as toroidal or helical magnetic field. Our 
simulations show that variable, low-density jets with velocities rsj 1000 km s 1 and r^j 10 times lighter 
than the environment can produce internal knots with significant synchrotron emission, and thermal 
X-rays in the shocked region of the leading bow shock moving in a dense medium. While models with 
a purely toroidal magnetic field show a very large degree of polarization, models with helical magnetic 
field show lower values and a decrease of the degree of polarization, in agreement with observations 
of protostellar jets. 

Subject headings: Herbig-Haro objects — ISM: jets and outflows — magnetohydrodynamics (MHD) 
— polarization — radiation mechanisms: non-thermal — shock waves 


1. INTRODUCTION 

Jets are present in astrophysical sources with various 
spatial scales, from Young Stellar Objects (YSOs) to ac¬ 
tive galactic nuclei (AGN). These collimated outflows are 
generally considered to be the result of bipolar ejection 
of plasma, associated with accretion onto a central ob¬ 
ject (Blandford & Payne 1982). Variability in the ejec¬ 
tion speed can produce internal shocks clearly seen in 
YSO jets in the form of bright optical knots (e.g., Raga 
& Noriega-Crespo 1998; Masciadri et al. 2002) called 
Herbig-Haro (HH) objects. 

Jets that arise from active galaxies can have relativistic 
speeds, and are well known synchrotron radiation emit¬ 
ters (see for example Tregillis et al. 2001; Laing et al. 
2006; Gomez et al. 2008). In contrast, YSO jets are non- 
relativistic and typically thermal radio sources. However, 
a few stellar sources, such as Serpens (Rodriguez et al. 
1989), HH 80-81 (Marti et al. 1995), Cepheus-A (Garay 
et al. 1996), W3(OH) (Wilner et al. 1999), and IRAS 
16547-4247 (Garay et al. 2003) present radio emission 
with negative spectral index interpreted as non-thermal 
(synchrotron) radiation. Notably, polarized radio emis¬ 
sion was detected in the jet of HH 80-81 (Carrasco- 
Gonzalez et al. 2010). Therefore, an interesting question 
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to answer is how jets with velocities of several hundreds 
km s -1 moving into a dense medium are able to produce 
shocks where particles can be accelerated up to relativis¬ 
tic energies and produce synchrotron radio emission. 

Synchrotron maps have been computed from MHD nu¬ 
merical simulations by several authors in different con¬ 
texts such as pulsar wind nebulae, e.g. Del Zanna et al. 
(2006), Volpi et al. (2008); supernova remnants, e.g. Or¬ 
lando et al. (2007); and accretion disks, e.g. Goldston 
et al. (2005). Broderick & McKinney (2010) and Porth 
et al. (2011) have performed MHD numerical simulations 
of relativistic jets. In particular, Porth et al. (2011) 
have performed MHD numerical simulation of relativistic 
AGN jets in order to study the synchrotron emission at 
the jet acceleration region by computing synthetic emis¬ 
sion maps of the spectral index, polarization degree and 
Rotation Measure (RM). In the case of non-relativistic 
(YSO) jets, given the large densities of such a jets at the 
launching region, the base of the jet is a thermal emitter. 
However, as the jets propagates the density decreases and 
non-thermal signatures can appear. 

We present a polarization study in order to shed light 
on the understanding of the non-thermal emission in pro¬ 
tostellar jets. We model, by using axisymmetric, mag¬ 
netohydrodynamic (MHD) simulations, the synchrotron 
emission, and we compute the resulting polarization 
map. The paper is organized as follows: in Section 2, 
we describe the model and the numerical setup; in Sec¬ 
tion 3 we show the results (synthetic radio, polarization, 
and X-ray emission maps); and in Section 4 we present 
our conclusions. 

2. NUMERICAL CALCULATIONS 
2.1. Initial setup 
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Fig. 1.— Coordinate system employed for simulating synchrotron emission. The rz-plane is the 2D plane of our axisymmetric simulation. 
The plane of the sky or image plane is the x'z '-plane and the y'-axis is the LoS. 


Our study is based on 2.5D axisymmetric, MHD sim¬ 
ulations carried out with the adaptive mesh refinement, 
eulerian code Mezcal (De Colle & Raga 2006; De Colle 
et al. 2008, 2012). We consider a 2D axisymmetrical 
adaptive grid, with a size of 0.2 and 0.5 pc along the 
r— and 2 —directions, respectively, and a maximum spa¬ 
tial resolution of 1.56 x 10 -4 pc, corresponding to 1280 
x 3200 cells (at the maximum resolution) along the r— 
and 2 —directions, and 6 levels of refinement. The en¬ 
vironment in which the jet propagates is homogeneous, 
with a uniform density n env = 3000 cm -3 , temperature 
T env = 100 K, and magnetic field Bq. At every timestep 
the jet is imposed in the computational domain by rewrit¬ 
ing the values of the physical parameters inside a region 
of the computational domain with r < Rj et = 0.03 pc 
and 2 < 0.003 pc, with density nj et = 300 cm -3 (Araudo 
& Rodriguez 2012) and velocity Uj et (along the 2 -axis). 
The longitudinal magnetic field (imposed on all compu¬ 
tational domain) is B z = B 0 , and the toroidal component 
is given by (Lind et al. 1989) 
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where R m = 0.02 pc, and B m is given in Table 1. In 
models M4 and M5 B z and B m are chosen so that B 
(= \JB 2 Z T L>m) results of the order of 0.1 mG (Carrasco- 
Gonzalez et al. 2010; Curran & Chrysostomou 2007). 
The jet pressure profile is constructed to ensure total 
pressure equilibrium at t = 0: 
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where /3 m — p env /(5 m /87r) and Penv — rz env fc^T 1 env . 

We consider five different initial configurations. Model 
Ml represents a continuous jet with constant injection 
velocity fj e t = no, whereas models M2-M5 have a time 
dependent injection velocity of the form 


Det = V 0 (l + Av cos (ut)), (3) 

where vq = 1000 kms -1 is the mean velocity of the 


flow, and uj = 27r/r, r = 50 yr and Av are the fre¬ 
quency, periodicity, and amplitude of the variability, re¬ 
spectively. The values of B z , Av and jet maximum ve¬ 
locity u max = uo(l T An) for the different models are 
given in Table 1. With these values, Uj et is in the range 
of ~ 600 — 1400 km s -1 , as observed in HH80-81 (Marti 
et al. 1995, 1998). 

We have also explored the case of a dense and slow jet 
(model M6). This model has the same parameters as the 
model M4, except that no = 300 km s -1 , An = 0.3, and 
the density of the jet and the surrounding medium are 
1000 and 100 cm -3 , respectively (see Table 1). 


2.2. Synthetic emission maps 
2.2.1. Non-thermal radio emission and Stokes parameters 

In this section, we present a brief description of the 
strategy used to compute the non-thermal synchrotron 
emission and the Stokes parameters. For details of the 
method we refer the reader to Ghisellini (2013). See also 
Rybicki & Lightman (1986) for an in-depth description 
of synchrotron emission. 

Synchrotron emission in YSO jets is produced by rela¬ 
tivistic electrons accelerated in internal and termination 
shocks (see Section 3.1). In the present study, we as¬ 
sume that there is a population of relativistic electrons 
with power-law energy distribution: 

n e = K 7 - p , (4) 

for 7 m i n < 7 e < 7max (being 7 e the Lorentz’s factor), 
and n e = 0 otherwise. We fix p = 2.1 in our calculations 
and determine K and 7 m i n assuming that the number 
density of non-thermal electrons is a fraction Xn (< 1) of 
n g , the electron density of the gas (assuming as well that 
the plasma is composed by the same number of electrons 
and protons, i.e. fully ionized hydrogen, in post-shocked 
regions): 

TTrnax — P+1 

Xn n g = / K^d le ~ K^-, (5) 

4 7m in P 4 


and that the energy density is a fraction (< 1) of the 
gas kinetic energy density e = m p n g v g/2: 
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Fig. 2. — Number density stratification maps, in units of 10 3 cm -3 , and displayed in linear color scale (see colorbar at the top). The 
black arrows depict the velocity field, with a scale shown at the bottom of the leftmost panel. The integration time is t = 1500 yr in all 
models. 
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Fig. 3. — Maps of the magnetic field intensity, (see colorbar at the top for the scale, in units of mG) obtained for all models at an 
integration time of 1500 yr. 
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TABLE 1 
Initial setup 


Model 

B z [mG] 

Bm [mG] 

^env/^jet 

Av 

%ax [km s - 1 ] 

Ml 

0 

0.1 

10 

0 

1000 

M2 

0 

0.1 

10 

0.2 

1200 

M3 

0 

0.1 

10 

0.4 

1400 

M4 

0 . 1/72 

0.1/72 

10 

0.4 

1400 

M5 

o.i/vTo 

o.3/7Io 

10 

0.4 

1400 

M6 

0.1/72 

0.1/72 

0.1 

0.3 

390 
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where m e is the rest mass of the electron and c is the 
speed of light. From equations (5) and (6), we obtain 

1 X e e 
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(8) 


We are interested in the study of the synchrotron radia¬ 
tion at frequencies v ~ 1 GHz. This emission is optically 
thin at frequencies larger than the self-absorption fre¬ 
quency (equation 4.56 of Ghisellini 2013) 


*4a = 


n 3 / 2 eRK „ , J 2/(p+4) 


4 B 


fa (p) 


integrated along the LoS (I{y) = f LoS j s (v)&y'), which 
was chosen to be y' (see Figure 1). In this study, a depen¬ 
dence with the angle between the shock normal and the 
post-shocked magnetic field is not considered (see for in¬ 
stance Orlando et al. 2007; Petruk et al. 2009; Schneiter 
et al. 2015 for a discussion of the acceleration mecha¬ 
nisms in supernova remnants). 

From the synchrotron emission we have also carried a 
polarization study by means of maps of the Stokes pa¬ 
rameters Qb and Ub , which can be computed as: 
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where z/l = eB/{2iim e c) is the Larmor frequency, e is the 
charge of the electron, and (equation 4.52 of Ghisellini 
2013): 
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/ a (p) -3 ^ ( + 4Q 
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( 10 ) 


By setting R = 10 17 cm as the size of the emitting region, 
we obtain z/ sa ~ 3 MHz 1 GHz for typical values of 
density and velocity in our simulations. 

In the optically thin case, the synchrotron specific in¬ 
tensity, for the isotropic case, can be written as (see equa¬ 
tion 4.45 Ghisellini 2013): 


(see e.g. Clarke et al. 1989; Jun & Norman 1996), where 
(x',z' are the coordinates in the plane of the sky (see 
Figure 1, dy' is measured along the LoS, <j)(y') is the 
position angle of the local magnetic field on the plane of 
the sky, and 


_ a + 1 
/0 “ a + 5/3 


(19) 


is the degree of linear polarization. The intensity of the 
linearly polarized emission is given by 


Ip{x', z', v) = yjQ%(x', z', v) + XJ%(x', z', v) (20) 
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where ot = 6.65 x 10 25 cm 2 is the Thomson cross sec¬ 
tion, up = B 2 /(8tt) is the magnetic energy density, and 


fj(p) ~ 3 p/2 (2.25 /p 2 - 2 + 0.105) (12) 

(see equation 4.46 Ghisellini 2013). Using equations (7) 
and (8), equation (11) gives: 


3,{v) = K sX 1 - 2a xl a n g v A g a B a + 1 u- a , (13) 

where B± is the component of the magnetic field perpen¬ 
dicular to the line of sight (LoS) and a = (p— 1)/2 = 0.55 
is the spectral index. The factor K s in equation (13) is: 


I<s = KiK 2 m\~ 3a c 2 ~ 5a (fim H ) 2a fj(p) (14) 


and the map of the position angle of the magnetic field 
(which gives the orientacion of the magnetic field in the 
plane of the sky) is computed as 

Xb(x',z') = ^tan~ 1 (U B {x',z , ,iy)/Q B (x',z',u)). (21) 

2.2.2. Thermal X-ray emission 

We calculated the thermal emission by integrating the 
free-free emissivity j u (n g ,T) along the LoS. In the low 
density regime, j u (n g ,T) = A(T), where n g is the 
electron density and T is the temperature. As with non- 
thermal radio emission, we assume that the post-shocked 
medium is fully ionized. The function A(T) was con¬ 
structed with the CHIANTI atomic database (Dere et al. 
1998) considering the energy range [0.15-8] keV and as¬ 
suming solar metallicity. 


where y is the molecular weight and ra# is the proton 
rest mass. The factors K\ and are 


K i 


ip~ !) 



(p- 1) 


(15) 


3. RESULTS 

The numerical simulations with the initial conditions 
summarized in Table 1 were carried out until they reach 
an integration time of 1500 yr. 


and 

K 2 = ^T^cr T (7re) (2a - 1)/2 . (16) 

Synthetic synchroton intensity maps are obtained from 
our 2D simulation in the following way. For each cell of 
our 2D axisymmetric grid we compute the synchrotron 
emissivity. The 2D plane of simulation (rz-plane) is 
tilted with respect to the plane of the sky (around the 
x'-axis) by an angle (p, and it is revolved around the 
symmetry axis (z), sampling a large number of angles in 
the 6 direction in order to obtain a “3D distribution” of 
the synchrotron emissivity. Then, the emissivity j s (v) is 


3.1. Shocks in proto stellar jets 

Figure 2 displays the number density stratification and 
the velocity field. As the jet interacts with the surround¬ 
ing medium, it forms a double shock structure, where the 
environment gas is accelerated by a forward shock, and 
the jet plasma is decelerated by a reverse shock. This 
structure, as well as the contact discontinuity separating 
the shocked interstellar material from the shocked jet 
material, is clearly visible at the head of the jet shown 
in Figure 2. In all cases, a slow bowshock travels against 
the surrounding environment with velocities ^ [200—260] 
km s -1 . 











6 


Cecere et al. 





o 

T—I 



123456 123456 123456 123456 123456 

x' [10 17 cm] x' [10 17 cm] x' [10 17 cm] x' [10 17 cm] x' [10 17 cm] 


Fig. 4.— 


Comparison of synthetic maps of the intensity of the linearly polarized radio emission obtained for all models at v — 5 GHz. 


Internal shocks are present only in the models where 
Av 7^ 0. Several jumps in axial velocity are present which 
mark the position of the internal shocks and have values 
in the range [400 — 500] km s -1 . With these values and 
the Rankine-Hugoniot jump conditions, we can estimate 
internal shocks velocities as large as 1000 km s _1 . These 
velocities, and considering that these shocks move in a 
low density medium, with densities of the order of 200 
cm -3 imply that the internal shocks are adiabatic 6 . Note 
that in models M3-M5 (Av = 0.4) the bow shock is 
significantly slower (200 km s -1 ) than in models Ml and 
M2 (260 km s _1 ). 

Figure 3 displays the magnitude of the magnetic field 
B. In models M2 and M3 the jet variability is quite evi¬ 
dent by the presence of a thin Mach disk in several inter¬ 
nal working surfaces. In contrast, the map corresponding 
to model M4 and M5 have a more complex morphology 
with less defined working surfaces, and a larger cocoon 
structure due to the magnetic field along the symmetry 
axis. 

3.2. Radio polarization 

Figure 4 shows synthetic maps of the intensity of the 
linearly polarized radio emission /p, at 5 GHz, for all 
models. The spatial resolution of these maps is 10 -3 pc. 
We have considered that the jet axis is tilted 15° with 
respect to the plane of the sky. 

In model Ml (Av = 0) the jet develops a single radio 
knot associated with material within the Mach disk (see 
Figure 4, left panel). The extended feature observed at 
z' ~ 6 x 10 17 cm, associated with an internal shock, is an 
artifact produced by the reflective boundary condition 
imposed along the symmetry axis. 

6 The cooling distance (see e.g. the equation (6) of Raga et al. 
2002) results larger than the jet radius, implying an adiabatic na¬ 
ture for the internal shocks. 


In contrast, the other models (Av ^ 0) display knotty 
radio structures produced by the internal shocks 7 . These 
knots decrease in brightness with distance from the jet 
source. Model M4 (B z = B m ) also shows radio emission 
in the region behind the main bow shock. The radio 
emission of this region in model M5, which has B z = 
B m /3 (see Table 1), results in a lower emission compared 
with model M4. 

In Figure 5 the degree of polarization of the syn¬ 
chrotron radiation (Ip(x f , z', v) / I(x ', z f ,v)) shows an im¬ 
portant result. Models M1-M3 exhibit a high degree 
of polarization of the synchrotron emission while model 
M4 displays strong variations toward the symmetry axis. 
This behavior is also observed in model M5, although to 
a lesser degree. These results can be understood consid¬ 
ering that, in an helical magnetic field, emission from re¬ 
gions with linearly polarized synchrotron radiation whose 
polarization directions are orthogonal to each other can¬ 
cel out when the emission is integrated along a LoS nearly 
perpendicular to the axis of symmetry of jet. A decrease 
of the degree of polarization has been observed in the 
jet associated with HH 80-81 (Carrasco-Gonzalez et al. 
2010 ). 

Figure 6 displays maps of the distribution of the posi¬ 
tion angle of the magnetic field xb(%',%')•> As in obser¬ 
vational studies, these xp(x / ,^ / )’s maps were performed 
from the synthetic maps of the Stokes parameters Qb 
and Ub, by using of equations (17), (18), and (21). These 
maps show that for models M4 and M5, which display 
variations in the degree of polarization, the orientation 
of the magnetic field in the plane of the sky has a compo- 

7 We are considering a jet with a fixed axis, and thus the working 
surfaces move along the axis of symmetry of the jet and never exit 
the cocoon carved by the main bow shock. Thus, their interaction 
is with the previously ejected jet material and not with the external 
medium. 
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Fig. 5.— Same as Figure 4 but displaying maps of degree of polarization 


nent parallel to the symmetry axis (given an additional 
support of the helical nature of the magnetic field), while 
in models M1-M3 the magnetic field is mostly perpendic¬ 
ular to it. 

A comparison between synchrotron and thermal X- 
ray emission for models M3 and M4 is shown in Fig¬ 
ure 7. For both models, the X-ray emission maps also 
display a knotty structure although less defined than 
its radio counterpart. Most of the thermal X-ray emis¬ 
sion comes from the environment material swept up 
by the main bow shock, as shown by Bonito et al. 
(2004, 2007, 2010) in hydrodynamic simulations. As it 
was mentioned above for the map of the linearly po¬ 
larized intensity, the total synchrotron emission maps 
display bright knots close to the central source. These 
bright knots emit a radio flux, at 5 GHz, of the order 
of 1.5 Xn~ 2a Xe a x 10 -18 erg s -1 sr -1 cm -2 Hz -1 , which is 
(aside of setting the exact value for the fractions Xn and 
Xe) in reasonable agreement with the flux reported by 
Carrasco-Gonzalez et al. (2010) for the knots in HH 80- 
81 (1 mJy/beam or 10 _18 erg s -1 sr -1 cm -2 Hz -1 ). 

Figure 8 shows that the synchrotron emission for the 
case of a dense and slow jet (model M6) is 30 times 
lower than the emission for a lighter and faster jet (model 
M4). Furthermore, it is important to do a comparison 
of the magnitude obtained for both the non-thermal and 
thermal emission mechanisms in radio wavelengths. The 
thermal radio-continuum emission is optically thin for 
the parameters chosen in our simulations 8 . Therefore, 
one can compute the ratio of non-thermal to thermal 
emissivities (using equation (5) of Velazquez et al. (2007) 
and equation (11) of this paper) in the bright radio knots 
of models M4 and M6 (located at z' ~ 1.8 x 10 17 cm and 
z' ~ 0.5 x 10 17 , respectively). For model M6 this ratio is 

8 The opacity r t h is <C 1, considering the equation (3) of 
Velazquez et al. (2007). 


js(^)/ith(^) = 0.03 at a frequency v = 5 GHz, while for 
model M4 the ratio is 1.4. Thus, a dense and slow jet does 
produce synchrotron emission, however at a level that is 
negligible compared to the thermal radio-continuum. 

4. DISCUSSION AND CONCLUSIONS 

Carrasco-Gonzalez et al. (2010) have shown the exis¬ 
tence of polarized radio emission associated with the HH 
80-81 protostellar jet. However, this issue has not been 
studied by HD or MHD simulations. We present the re¬ 
sults obtained from 2.5D MHD simulations of low- and 
high-density YSO jets. We have considered cases with 
a constant, and a time-dependent jet ejection velocity. 
Furthermore, cases in which the magnetic field is toroidal 
or helical were analyzed. 

Assuming a population of relativistic electrons which 
are accelerated within stellar shock jets, we have used 
standard prescriptions to estimate their synchrotron 
emission. Our results indicate that while the thermal X- 
ray emission is dominated by the shocked environment 
material at the head of the jet (in agreement with Raga 
et al. 2002; Bonito et al. 2004, 2007, 2010), the non- 
thermal radio emission turns out to be dominated by the 
jet material inside the internal shocks. 

Also, radio maps reveal that the variability in jet veloc¬ 
ity is important in generating bright knots of synchrotron 
emission, produced when slow jet material is caught up 
by faster jet material. 

Our models show that a jet with a toroidal magnetic 
field emits synchrotron radiation with a high degree of 
polarization. In contrast, models with a helical magnetic 
field exhibit a decrease on the degree of polarization, 
in good agreement with observational results (Carrasco- 
Gonzalez et al. 2010). 

Finally, our results indicate that non-negligible syn¬ 
chrotron emission can be obtained in low-density and 
high-velocity protostellar jets. 


















Cecere et al. 

01 23\ 45\ 68^ 90- [°] 


o 


Ml 

M2 

M3 

M4 

M5 

p 

r 


in** 



; 

n 

m 

w \ 


r 

ms ) 

-** 

ft' * 


i 

i 

1 

s 

{ 


w 

J 

1 

s 


¥ 




\ 

* 


* 

V 

\ 

m 

m 



1 2 3 4 5 6 

x' [10 17 cm] 

1 2 3 4 5 6 

x' [10 17 cm] 

1 2 3 4 5 6 

x' [10 17 cm] 

1 2 3 4 5 6 

x' [10 17 cm] 

1 2 3 4 5 6 

x' [10 17 cm] 


Fig. 6. — Same as Figure 4 but showing the position angle of the magnetic field B (which is measured with respect to the z' — axis, as 
indicated by the stick marks in the colorbar). 
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Fig. 8.— Comparison of synthetic synchroton emission maps, at 5 GHz, of model M4 (left panel) with model M6 (right panel). The 
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